/*
 *  commonutil.h
 *  EPPI-0.0
 *
 *  Created by Juan Carlos Garcia Orden on 15/11/07.
 *
 */

#ifndef _commonutil_h_
#define _commonutil_h_

#include <string>
#include <vector>
#include <math.h>
#include <iostream>
#include "MATH/sparsematrixslu.h"
#include "MATH/TENSOR/tensor.h"
#include "MATH/TENSOR/xtensor.h"
#include "MATH/TENSOR/tutils.h"
#include "Eigen/Dense"

using namespace std;

	void addSystemName(string& str, const string& strname);

	string getStringHead(const string& str);

	string getStringTail(const string& str);

	string getStringExtension(const string& str);

	bool checkForSlash(const string& str);

	string subsSlash(const string& str);

	double lTwoNorm(const vector<double>& v1);
	
	double lTwoNorm(const vector<double>& v1, const int&);
	
	double lTwoNorm(const vector<double>& v1, const vector<double>& v2);

	double lTwoNormSq(const vector<double>& v1);

	void printSTLVector(vector<double>& v);
	
	void printEigenVectorXd(const VectorXd & v);

	void printSTLVector(const vector<string>& v);
	
//	void printSTLVector(const vector<int>& v);
	
	void printSTLVector(vector<int>& v);
	
	void uVersor(const vector<double>& v1, const vector<double>& v2, const double & dis, vector<double>& u);
	
	void uVersor(const vector<double>& v1, const double & dis, const int& dim, vector<double>& u);	

	void printSTLMatrix(const vector< vector<double> >& m);

	void printSTLMatrix(const vector< vector<int> >& m);
	
	void printEigenMatrixXd( const MatrixXd & m);
	
	void printEigenMatrix3d( const Matrix3d & m);

	void printSparseMatrixSLU(sparseMatrixSLU& m);

	void printSparseMatrix(sparseMatrix& m);

	void setZero(vector<double>& v);

	void setZero(vector< vector<double> >& m);
	
	void lu(vector<vector<double> >& mat);	// Replaces matrix mat by its LU decomposition
//Arma	void lu(mat &);	// Replaces matrix mat by its LU decomposition
	
	void forSub(vector<vector<double> >& mat, vector<double>& vec);	// Perform a forward subs, mat already factorized. Vec is replaced
//Arma	void forSub(mat & mat, rowvec & vec);	// Perform a forward subs, mat already factorized. Vec is replaced

	void forSub(vector<vector<double> >& mat, vector<vector<double> >& mat_nrhs);	// Perform a forward subs, mat already factorized. mat_nrhs is replaced
	void forSubTrans(vector<vector<double> >& mat, vector<vector<double> >& mat_nrhs);
//Arma	void forSubTrans(mat &, mat &);
	
	void backSub(vector<vector<double> >& mat, vector<double>& vec);	// Perform a backward subs, mat already factorized. Vec is replaced
//Arma	void backSub(mat &, rowvec &);	// Perform a backward subs, mat already factorized. Vec is replaced

	void backSub(vector<vector<double> >& mat, vector<vector<double> >& mat_nrhs);	// Perform a backward subs, mat already factorized. Vec is replaced
	void backSubTrans(vector<vector<double> >& mat, vector<vector<double> >& mat_nrhs);
//Arma	void backSubTrans(mat &, mat &);

	int locateInSparse(vector<int>& ia, vector<int>& ja, int& I, int& J);

	void inverse(vector<vector<double> >& mat);	// Inverse based on a LU decomposition

	void rigidTransOver(vector<double>&, vector<double>&, vector<vector<double> >&);	// Rigid transformation. x is overwritten

	void rigidTrans(vector<double>&, vector<double>&, vector<vector<double> >&, vector<double>&);

	void getGlobalCoords(vector<double>&, vector<double>&, vector<double>&);

	void getLocalCoords(vector<double>&, vector<double>&, vector<double>&);

	void getGlobalCoordsNoTrans(vector<double>&, vector<double>&, vector<double>&);

	void substractSTLVectors(vector<double>& a, vector<double>& b);	 // a-b, a is overwritten

	void addSTLVectors(vector<double>& a, vector<double>& b);	 // a+b, a is overwritten

	void crossProdSTLVectors(vector<double>& a, vector<double>& b);  // a x b, b is overwritten

	void evalInBeta(double& beta, vector<double>& vBeta, vector<double>& vN, vector<double>& vNPlus);

	void evalMid(vector<double>& vMid, vector<double>& vN, vector<double>& vNPlus);
	
	void VoigtFormat(const vector< vector < double > >& T2, vector < double > & T2v);
	
	void VoigtFormat(const Tensor2 & T2, vector < double > & T2v);
	
	void VoigtFormat(const Matrix3d & T2, VectorXd & T2v);
	
	void VoigtFormat(const Tensor2 & T2, VectorXd & T2v);
	
	void VoigtFormat(const vector< vector < vector < vector < double > > > >& T4, vector < vector < double > > & T4v);
	
	void VoigtFormat(const Tensor4 & T2, vector < vector < double > > & T4v);
	
	void VoigtFormat(const Tensor4 & T2, MatrixXd & T4v);
	
	void InverseVoigtFormat(Matrix3d & Cf, const VectorXd & C);
	
	void InverseVoigtFormat(Tensor2Sym & CT, const VectorXd & C);
	
	double dotVoigt( const VectorXd & V1, const VectorXd & V2);
	
//	VectorXd dotVoigt( const MatrixXd & M1, const MatrixXd & M2);
        
	VectorXd makeVoigt(const VectorXd & C);

    VectorXd makeCiVoigt(const VectorXd & C);

	SuperMatrix convert_stl_SLU(vector<vector<double> >&);

	void convert_SLU_stl(SuperMatrix&, vector<vector<double> >&);

#endif
